Cite this Article
Zhang Shujun, Wang Shaofeng. Modification of the Peierls–Nabarro model for misfit dislocation. Chinese Physics B, 2020, 29(5): 056102  
Permissions
Modification of the Peierls–Nabarro model for misfit dislocation
Zhang Shujun, Wang Shaofeng †
Department of Physics and Institute for Structure and Function, Chongqing University, Chongqing 400030, China

 

† Corresponding author. E-mail: sfwang@cqu.edu.cn

Project supported by the National Natural Science Foundation of China (Grant No. 11874093).

Abstract

For a misfit dislocation, the balance equations satisfied by the displacement fields are modified, and an extra term proportional to the second-order derivative appears in the resulting misfit equation compared with the equation derived by Yao et al. This second-order derivative describes the lattice discreteness effect that arises from the surface effect. The core structure of a misfit dislocation and the change in interfacial spacing that it induces are investigated theoretically in the framework of an improved Peierls–Nabarro equation in which the effect of discreteness is fully taken into account. As an application, the structure of the misfit dislocation for a honeycomb structure in a two-dimensional heterostructure is presented.

PACS: ;61.72.Bb;;61.80.-x;
Keyword:;interfacial misfit dislocation;;the energy of misfit dislocation;
1. Introduction

Dislocations are structural defects that are central to the understanding of the mechanical properties of materials. With the recent development of optoelectronic devices based on heterostructures, it is crucial to understand the properties of the misfit dislocations that occur at the interface between two materials with different lattice constants. In fact, this is a topic that has received much attention for more than fifty years. To present the basic physical concept of misfit dislocation in the simplest way, Frank and Van der Merwe[1] extended earlier work on the Frenkel–Kontorova (F–K) model,[2] which describes an atomic chain linked by nearest-neighbor harmonic forces and subjected to a periodic substrate potential, taking account of the difference in lattice spacing between the chain and substrate. However, the dislocation obtained by the F–K model exhibits a short-range deformation similar to an exponential decay. On the basis of continuum elasticity theory, Dundurs[3] obtained the Airy stress function for an edge dislocation, using an analogy between concentrated forces and edge dislocations.[4] As in the original Peierls–Nabarro (P–N) model,[5,6] the misfit dislocations result from a nonlinear interaction between two semi-infinite crystals. Based on an extension of the P–N theory by van der Merwe,[7] Yao et al.[8] proposed an equation for an interfacial misfit dislocation that was similar to the P–N equation. By using an exact solution given first by van der Merwe[9] and later recovered by Yao et al.,[8] it is possible to use continuum elasticity theory to describe the long-range deformation of a dislocation at length scales outside the core field of the dislocation line. However, in general, it is clearly unrealistic to view a solid as a continuous medium, ignoring the discrete nature of the crystal lattice. Because the dislocation core structures play an important role in many crystal plasticity phenomena, it is of great interest to describe these structures accurately on the atomic scale.

Recent theoretical analyses of dislocation configurations have been based on molecular dynamics simulations[10] or local density function simulations. It is possible to obtain a misfit P–N equation by solution of the balance problem for a semi-infinite lattice,[11] and the resulting equation is similar to that derived by Yao et al.[8] However, in this approach, following that in the original P–N model, when the upper and lower semicrystals are glued together along their interface, the effect of the lower semicrystal on the bottom surface of the upper semicrystal is ignored, as is the effect of the upper semicrystal on the top surface of the lower semicrystal. Furthermore, the changes in the interfacial spacing introduced by misfit dislocations have also been ignored in previous studies.

In this paper, a correction to the misfit dislocation equation that takes account of discreteness is derived by modifying the displacement fields that satisfy the balance equations. An extra term proportional to the second-order derivative, which originates from a discrete feature, it is used to reflect the deformation that occurs near to the dislocation core. It is reasonable to view the modified equation as a combination of the F–K equation and the P–N equation. The integral term describes the elastic effect of the interior of the semi-infinite lattices regarded as continuous media, while the differential term describes the short-range interaction between atoms at the misfit interface and represents the effect of discreteness. To estimate the change in interfacial spacing induced by misfit dislocation, we consider the variation in the generalized stacking fault energy along the direction of interfacial spacing.

2. Discreteness correction to misfit dislocation equation

As in the original P–N model, misfit dislocations are created by gluing together two semi-infinite crystals of different materials, with each part being a regular crystal with the cut plane as its surface. As shown in a previous investigation of the balance problem for a semicrystal acted on by external forces on the surface,[12] the equilibrium equation for the atoms on the surface is

where the displacement uj and the force per unit area fi (i,j=x,y) are two-dimensional (2D) vectors, Λ is a 2 ×2 matrix, and we have used the Einstein convention, whereby repeated indices in the same expression are summed over and the summation sign is omitted. We are then able to obtain the misfit P–N equation in continuum elasticity theory.[11] This equation describes the long-range deformation of the dislocation for length scales outside the core field of the dislocation line. However, it is not applicable to the deformation at an atomic scale in the neighborhood of the dislocation core.

In the continuum approximation, it is difficult to distinguish between the surface and the interior of a semi-infinite crystal, and the surface effect is ignored. The second-order derivative describes the lattice discreteness effect that arises from the surface effect. To represent the effect of the discrete crystal lattice and to take account of the short-range interaction between the interfacial atoms, we consider a correction based on a reduced dynamical matrix for the surface of a semi-infinite lattice.[13] Taking into account that a heterostructure is composed of two dissimilar crystals, it is necessary to distinguish between the lattice parameters of materials 1 (upper) and 2 (lower), denoted by a and b, respectively, with . The shear moduli of the upper and lower semi-infinite materials are denoted by μa and μb, respectively, and their Poisson’s ratios by νa and νb. Then, the bottom surface of the upper semicrystal is termed the a plane, and the top surface of the lower semicrystal the b plane. The displacement fields of the planes are denoted by ua(x′) and ub(x′). We introduce the relative displacement and mass center displacement

On the basis of the rotational transformation between the respective coordinate systems for the upper and lower semicrystals, and by choosing the a plane as the unified direction of the common coordinate system, the balance equations for the modified displacement fields are transformed to wave-vector space (k space) as

with

where

The parameters βa and βb representing corrections for the effects of discreteness depend on the arrangement of atoms and are related to variations in bond length and angle, details of which may be found in Refs. [13,14]. The above matrices include two sets of terms: (i) those terms including k derived using elastic continuum theory and arising from long-range interactions between atoms in the interior of the semi-infinite lattice; (ii) terms depending on k2 that arise from the short-range interactions of atoms in the misfit interface and represent corrections for the effects of discreteness. From Eqs. (2) and (3), the relationship between the forces and the relative displacements is obtained as

where the matrix is given by

The Taylor series expansion of this matrix in the long-wave approximation with higher-order terms being neglected gives

The expansion coefficients are constants that depend on the properties of the two materials and are given by

with

The modulus operators in k space are

Fig. 1. Misfit dislocation in a BN/AlN heterostructure.
Fig. 2. (a) Energy per unit area, γ-potential (d0=2.171 Å), as a function of the relative slip and the change in interfacial spacing in the case of δx=0.6503 and δy=-4.5414. (b) γ-potential as a function of . (c) γ-potential as a function of the change in interfacial spacing . The solid curves are given by Eq. (22) and the discrete points are the results of a first-principles calculation.
Fig. 3. (a) Comparison of the left-hand side (red dashed curve) and right-hand side (blue solid curve) of Eq. (9) for fitting parameters g = 1.3735, c1=0.0569, and c2=0.0985. (b) Application of the dislocation density to compare the differences between the solution of the equation of Zhang et al.[11] and that of Eq. (9). The blue dashed curve represents the dislocation density of Eq. (9) in the case  eV, and δx=0, and the solid curve represents the dislocation density from Eq. (9) in the case  eV, and δx=0.6503. (c) Comparison of the left-hand side (red dashed curve) and right-hand side (blue solid curve) of Eq. (10) for fitting parameters cy1=0.040, cy2=-0.3833, and g1=1.74. (d) Change in interfacial spacing induced by misfit dislocation.

Equation (4) in k space can then be transformed into the following equations in real space:

Taking into account the relationship between the plastic displacement and the misfit displacement, the relative displacement between corresponding atoms on each side of the interface, , can be written as

where (as in the original P–N theory) is the displacement vector in the glide plane of the reference lattice. Taking account of the fact that for a periodic misfit dislocation both are periodic functions and using the identity a,[15] equations (7) and (8) can be simplified as follows:

In comparison with the equation of Yao et al.[8] and Zhang et al.,[11] Equation (9) includes the second-order derivative corrections to the result of the continuum elastic theory. It should be noted that these corrections provide a good description of the deformation around the dislocation core compared with the continuum result. Equation (10) determines the effect of the change in interfacial spacing when periodic misfit dislocations exist at the interface.

Using Eq. (2) together with Eq. (3), we can express the displacement of the mass center in terms of the relative displacement

where the matrix is given by

Similarly to , we can expand in a power series and neglect the higher-order terms,

where the expansion coefficients are given by

In view of the forms of the relative displacement in real space, equation (11) can be transformed into the following equations in real space:

with

Taking account of the fact that the relative displacements of interfacial atoms are periodic functions and again using the identity a,[15] we can simplify Eqs. (14) and (15) as follows:

3. Application

We consider a heterostructure in which the misfit occurs in a single direction, as illustrated by the 2D honeycomb structure in Fig. 1. As pointed out by van der Merwe,[9] a reference lattice with lattice parameter c is generated from lattice parameters a and b, and is defined by

where P is a positive integer. Thus, c and p (the period between dislocations) are given by

For a periodic misfit dislocation in which the upper and lower semi-infinite crystals are boron nitride (BN) and aluminum nitride (AlN) in a 2D honeycomb structure, the lattice parameters are chosen to be a = 2.513 Å and b = 3.127 Å, Poisson’s ratios νa=0.229 and νb=0.461, and shear moduli (as obtained from a first-principles calculation) μa=111.44  N/m and μb=38.92 N/m. For the reference lattice, Equation (19) then gives the values of the parameter c = 2.789 Å and the period p = 12.55 Å. According to the arrangement of the atoms in the heterostructure, the discreteness correction parameters are[14]

From the lattice parameters and elastic constants, we can then obtain the parameters in Eqs. (9) and (10), as shown in Table 1.

Table 1.

Parameters (N/m) in Eqs.(9) and (10).

.

In fact, the structure and properties of the dislocations depend on the restoring forces defined as the gradients of the generalized stacking fault energy.[16] Recent first-principles calculations of the generalized stacking fault energy (also generally called the γ-surface) have provided reliable values for the restoring forces. Because the lattice parameters of the materials do not match, as in the original P–N theory, we still assume that the effective interaction between two dissimilar semi-infinite crystals can be approximated using a reference lattice. The adhesion energy of the reference lattice varies with the interfacial distance, and this relationship can be approximated well by[11]

where e0 is the formation energy for the reference lattice when y=d0. The parameters ξ, Vs, and Vl are determined by the results of numerical fitting. To calculate the γ-surface for the heterostructure, we find a reference lattice with parameter c and distance d0 between layers in the glide direction. Taking into consideration the change of spacing in the interface caused by the periodic misfit dislocations, the γ-surface is generalized to a γ-potential that includes the changes in the interaction energy as the spacing varies.[12] The γ-potential can be regarded as the interaction energy between the two semi-infinite crystals. By definition, when the reference lattice is cut into two semi-infinite crystals by a given cut plane, the change in interaction energy comes from a relatively rigid translation between the two crystals. The rigid translation can be parallel or perpendicular to the cut plane, and so the interaction energy per unit area as a function of the translation is defined as the γ-potential. According to the relationship between the interfacial distance and energy, the γ-potential can be expressed as

where δx and δy are dimensionless parameters that are used to correct the generalized stacking-fault energy along the cut plane and interfacial spacing, respectively. The change in interfacial spacing caused by the periodic misfit dislocation is relatively small, therefore it is usual to expand the γ-potential in a Taylor series of this change in interfacial spacing as follows:

where the energy of adhesion .

To estimate the size and shape of the change in interfacial spacing caused by the misfit dislocation between the misfit planes, it is necessary to calculate the γ-potential of the reference lattice. These calculations are carried out with the Vienna ab initio simulation package (VASP), which is based on Kohn–Sham density functional theory. The generalized gradient approximation (GGA) as proposed by Perdew, Burke, and Ernzerhof (PBE) for the exchange–correlation functions[17] is used and the core electrons are dealt with using the projector augmented wave (PAW) approach.[18,19] The plane-wave cutoff is 600 eV in all the calculations. The criterion to halt the relaxation of the electronic degrees of freedom is set to be when the total energy change is smaller than 10−6 eV. The convergence iteration is halted when the forces on the ions are smaller than . The k points are chosen as 17 ×1 ×1 for the reference lattice. For BN/AlN, the k points are chosen as 3 ×1 ×1 Monkhorst–Pack meshes owing to the large size of the supercells. There is a 20 Å thick vacuum region to reduce interactions between inter-layers in the single-layer case. It should be noted in particular that we add a hydrogen bond saturation treatment to eliminate the boundary effects.

For an interfacial distance of the reference lattice d0=2.171 Å and an energy of adhesion  eV, figure 2(a) shows the results for the γ-potential (energy per unit area) calculated using Eq. (22). The corresponding parameters δx andδy are determined by fitting to the numerical results, as shown in Figs 2(b) and 2(c).For a given γ-potential, the restoring force per unit area can then be obtained as

The lattice parameter of the upper semi-infinite crystal is smaller than that of the lower semi-infinite crystal, therefore there exists a compressive stress at the interface. The dislocation is characterized by a slip distribution that satisfies the boundary conditions . It is difficult to obtain an exact solution of the improved P–N equation with full consideration of the effect of discreteness. However, referring to the solution of the original P–N equation and considering that the change in interfacial spacing has virtually no effect on the solution of Eq. (9), we can be obtain the following solution:[20]

with

where c1, c2, and g are parameters that need to be determined. On substituting this solution into Eq. (9), the left-hand side of the equation can be written as a polynomial function of y with q as its prefactor. On substituting Eq. (23) into Eq. (9), a similar result can be obtained for the right-hand side by a series expansion in c1 and c2 as small corrections. Taking into account that is small and can be ignored, we obtain the values of the parameters c1, c2, and g by truncating the series rationally using orthogonal polynomials,[21,22] as shown in Fig. 3(a).

In the improved P–N equation, the second-derivative term describes the lattice discreteness effect that results from the surface effect. To facilitate comparison of the differences (see table 2) between the solution of the equation of Zhang et al.[11] and that of Eq. (9), which includes corrections to the relationship between the relative displacement of the adhesive surface and the restoring stress, we define (as in the original P–N theory) the dislocation density ρs as the local gradient of the slip distribution ,

In Fig. 3(b), we use this dislocation density to illustrate the differences between the solutions of the two equations. The core structure of a dislocation predicted from the improved equation [Eq. (9)] is different from that predicted from the equation of Zhang et al. The characteristic length scale of misfit dislocation is defined as

For the solution of the equation of Zhang et al.,[11] the characteristic scale ζ=3.99 Å when the parameter g = 1.8306, whereas for the solution of Eq. (9), the characteristic scale ζ=5.03 Å when g = 1.3735. This clearly illustrates the considerable increase in the characteristic scale ζ of misfit dislocation when the discreteness corrections are included.

Table 2.

Parameters g, c1, and c2 in the misfit dislocation solution.

.

We can estimate the change in interfacial spacing on the basis of Eq. (10) with the given dislocation solution . Here we assume the following form for the change in interfacial spacing induced by the misfit dislocation:

with

The dimensionless parameters g1, cy1, and cy2 can be determined by fitting to Eq. (10), as shown in Fig. 3(c). These parameters are used to obtain the relative displacement induced by the misfit dislocation along the direction of interfacial spacing, which is shown in Fig. 3(d).

From the lattice parameters and elastic constants, we obtain the dimensionless expansion coefficients in Eqs. (16) and (17), as shown in table 3. The interfacial displacements can be obtained from the relationship between the displacement of the mass center and the relative displacement. Figures 4(a) and 4(b) clearly show the differences in the interfacial displacement before and after the discreteness correction to the misfit dislocation equation. It can be seen that the theoretical result is closer to the numerical result for the displacement along the direction of the slip plane. Although the variation in interfacial spacing is taken into account when there are periodic misfit dislocations at the interface, the numerical result is about four times the theoretical prediction for the displacement along the direction of interfacial spacing. These discrepancies may originate from the approximation of the effective interaction when we take the equilibrium distance of the interface as d0 in the BN/AlN heterostructure. In Fig. 4(c), we compare the theoretical predictions for the deformation of the interface caused by misfit dislocation with the results of first-principles simulations after full relaxation.

Fig. 4. (a) Interface displacements along the x direction: the green and blue solid curves represent the theoretical predictions for taking account of the discreteness correction to the misfit dislocation equation, the green and blue dashed curves represent the theoretical results without these corrections, and the green and blue filled circles represent the respective results of the numerical simulation. (b) Comparison between theoretical and numerical results for the interface displacement along the y direction: the green and blue solid curves represent the theoretical predictions for , and the green and blue filled circles represent the numerical results. (c) Comparison between theoretical and numerical results for the interfacial atomic configuration: the green and blue open circles represent the theoretical predictions, and the red filled circles represent the results of a numerical simulation based on a first-principles calculation after full relaxation.
Table 3.

Dimensionless coefficients in Eqs. (16) and (17).

.

In fact, the equilibrium slip distribution is determined by the competition between two contributions to the energy of misfit dislocation. The first is the elastic strain energy stored in each semicrystal, which can be obtained through the variational formulation of the P–N theory by adding up the contributions from each infinitesimal dislocation and can be written as

with

where . The first term EL in the elastic strain energy comes from the long-range interaction and can be determined using elastic continuum theory. The second term Es comes from the contact interaction resulting from the discreteness effects at the interface. The second contribution to the energy of misfit dislocation takes account of the atomic interaction across the glide plane and reflects the restoring stress that is balanced by the sum of the contributions from the infinitesimal dislocations. It can be written as

To determine the energy of the misfit dislocation, we introduce the adhesion work, which can be written as

An analogous definition may be given for the misfit interface based on first-principles calculations, in which case the heterostructure must be divided into two semi-infinite crystals along the interface to determine the adhesion work, which can then be calculated as

where is obtained from the first-principles calculation after full relaxation, and is the energy required for the rigid separation of the heterostructure at the interface in the absence of interaction between the two semi-infinite crystals when the distance of the interface is chosen as d = 6 Å. In fact, the adhesive work is one of the basic characteristics of interfacial bonding and is crucial to understanding the microscopic mechanism controlling the strength and fracture of composite materials. Misfit dislocations occur at interfaces in composite materials with lattice mismatch, and the associated energy can determine the adhesive work for the materials. As can be seen from table 4, correction for the effects of discreteness reduces the calculated value of the misfit dislocation energy. For the BN/AlN heterostructure, the adhesive work according to the first-principles calculation is -1.0569 eV/Å, compared with the theoretical value of -1.1690 eV/Å.

Table 4.

Energy of misfit dislocation (eV/Å ).

.
4. Summary

In an extension of the original P–N model, we have considered the effects at the interface when the upper and lower semicrystals are glued together. In general, these effects can be represented by second-order derivative terms to the F–K model. These terms can be viewed as corrections to the misfit dislocation equation to take account of the short-range interaction between atoms at the interface. The corrections presented here allow more accurate theoretical predictions of the interfacial atomic configuration. Using a γ-potential that includes the interaction energy due to the variation in interfacial spacing, together with an equation for the normal formation of the glide plane, we have estimated the change in spacing for an interface exhibiting misfit dislocation. We have presented an application of the improved P–N equation to an interfacial misfit dislocation in a BN/AlN heterostructure. Based on an approximate solution of the misfit dislocation equation, the theoretical result for the interface displacement along the glide direction is close to the numerical result. Although we have taken into account the change in interfacial spacing induced by misfit dislocation, it can be seen that there are some discrepancies between the theoretical and numerical results for the interface displacement along the y direction, and so the agreement between the results of analytical theory and those of numerical calculation needs to be further improved. We have found that the correction for the effects of discreteness reduces the value of the energy of misfit dislocation.

Reference
[1] Frank F C Van der Merwe J H 1949 Proc. R. Soc. Lond. A 198 205
[2] Frenkel Y I Kontorova T 1938 Zh. Eksp. Teor. Fiz. 8 1340
[3] Dundurs J Hetenyi M 1961 J. Appl. Mech. 28 103
[4] Dundurs J 1968 J. Appl. Phys. 39 4152
[5] Peierls R 1940 Proc. Phys. Soc. 52 34
[6] Nabarro F R N 1947 Proc. Phys. Soc. 59 256
[7] Van der Merwe Q 1963 J. Appl. Phys. 34 117
[8] Yao Y Wang T Wang C 1999 Phys. Rev. B 59 8232
[9] Van der Merwe J H 1950 Proc. Phys. Soc. A 63 616
[10] Yu T Xie H X Wang C Y 2012 Chin. Phys. B 21 026104
[11] Zhang S J Wang S F Wang C Y 2020 J. Appl. Phys. 127 085303
[12] Wang S F 2015 Philos. Mag. 95 3768
[13] Wang S F 2009 J. Phys. A: Math. Theor. 42 025208
[14] Zhang H L Wang S F Wang R Jiao J 2010 J. B Eur. Phys. 73 489
[15] Hirth J P Lothe J 1982 Theory of Dislocations (Krieger
[16] Rodney D Ventelon L Clouet E Pizzagalli L Willaime F 2017 Acta Mater. 124 633
[17] Perdew J P Burke K Ernzerhof M 1996 Phys. Rev. Lett. 77 3865
[18] Blöchl P E 1994 Phys. Rev. B 50 17953
[19] Kresse G Joubert D 1999 Phys. Rev. B 59 1758
[20] Wang S F Zhang S J Bai J H Yao Y 2015 J. Appl. Phys. 118 244903
[21] Wang S F 2003 Phys. Lett. A 313 408
[22] Wang S F Li S R Wang R 2011 J. B Eur. Phys. 83 15